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Abstract 

In  this  paper  we  establish  a  Joint  central  limit  theorem  for  customer  and  time  averages  by 
applying  a  martingale  central  limit  theorem  in  a  Markov  ftamework.  The  limiting  values  of  the 
two  averages  appear  in  the  translation  terms.  This  central  limit  theorem  helps  to  construa 
confidence  intervals  for  estimators  and  perform  statistical  tests.  It  thus  helps  determine  which 
finite  average  is  a  mote  asymptotically  efficient  estimator  of  its  limit  As  a  basis  for  testing  for 
PASTA  (Poisson  arrivals  see  time  averages),  we  determine  the  variance  constant  associated  with 
the  central  limit  theorem  for  the  difference  between  the  two  averages  when  PASTA  holds. 

OR/MS  Subject  ChMriflcatioa:  (Queues,  limit  theorems:  central  limit  theorems  for  customer 
and  time  averages.  Queues,  statistical  inference:  estimating  customer  and  time  averages. 
Statistics,  estimation:  averages  over  time  and  at  embedded  points. 
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This  paper  is  concerned  with  estimating  limiting  time  averages,  customer  (embedded) 
averages  and  their  differences.  When  the  two  limiting  averages  are  known  to  agree,  i.e.,  when 
arrivals  see  time  averages  (ASTA);  see  Melamed  and  Whitt  (1990a)  and  references  cited  there,  we 
want  to  know  which  finite  average  is  a  more  asymptotically  efficient  estimator  (i.e.,  produces 
snudler  confidence  intervals  with  large  samples).  When  the  two  limiting  averages  need  not  agree, 
we  want  to  estimate  their  difference  and  be  able  to  test  for  ASTA. 

To  illustrate  what  happens  when  ASTA  is  known  to  hold,  we  give  two  examples. 

Example  1:  The  Workload  In  the  M/M/1  Queue. 

Let  { f/(()  :  r  >  0}  be  the  continuous-time  workload  (or  virtual  waiting  time)  process  in  an 
M/M/1  queue  with  service  rate  I  and  arrival  rate  p  <  1.  Let  ( N(r)  :  r  >  0 )  be  the  Poisson  arrival 
process  with  associated  arrival  Umes  j  T,  ;  n  >  1 }.  Then  { (/(T,  -)  :  n  S  1 }  is  the  sequence  of 
waiting  times  (before  beginning  service).  It  is  well  known  that  the  time  average 

V(f)  =  -5- J' f/(5)ds.  f>0.  (1) 

t  '0 

and  the  customer  average 

1 

I  .  ^(0  >  0 .  (2) 

^(0  tTi 

both  convttge  with  probability  one  ( w.p.  1 )  to  p/(  1  -  p )  as  r  ee,  so  that  we  have  ASTA. 


□  □ 
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/ -w)  =*►  N{0,ai) 


(4) 


as  r  -» eo,  where  v  =  w  =  p/(l  -  p),  denotes  convergence  in  distribution  and  N(m,  <t~) 
denotes  a  normally  distributed  random  variable  with  mean  m  and  variance  o‘.  Moreover, 


o 


2 

V 


2p(3  -p) 
(1  -  P)^ 


and 


(2  +  Sp-4p‘  +  p^) 
(1  -P)^ 


see  Table  2  of  Whitt  ( 1989),  so  that 


(5) 


Cm 


OzML^^o. 

(I  -p)^ 


(6) 


Hence,  for  this  example,  the  time  average  V(t)  is  always  asymptotically  more  efficient.  (This 
result  may  seem  to  contradict  results  for  the  M/G/I  model  by  Law  ( I97S);  he  showed  tha.  W{t)  is 
superior  to  several  indirect  estimators,  but  V(r)  was  not  an  alternative  considered  there.) 

Example  2:  A  Uniformized  Continuous«Tiine  Markov  Chain. 


Let  {X(r)  :  r^O)  be  an  arbitrary  irreducible  finite-state  continuous-time  Markov  chain, 
uniformized  so  that  the  time  between  successive  transitions  has  an  exponential  distribution 
iixlependent  of  the  state.  (From  the  general  case,  this  is  achieved  by  introducing  extra  transitions 
from  states  to  themselves;  e.g.,  see  Keilson  (1979).)  Let  [N{t)  :  /  ^  0}  be  the  Poisson  process 
counting  the  transition  epochs  T„.  Then  {X(rn-}  :  n  ^  1 }  is  the  embedded  discrete-time 
Markov  chain.  Let  /  be  an  arbitrary  function  on  the  state  space  of  the  chain  and  let 
U{t)  =  /tX(01.  It  is  well  known  that  both  V{t)  in  (1)  and  W{i)  in  (2),  defined  in  terms  of  UU) 
and  T„  here,  converge  w.p.l  as  r  eo  to  £/[X(eo)]  where  X(ee)  has  the  unique  stationary 
distribution  of  {X(/)  :  /  2  0},  Moreover.  (3)  and  (4)  hold  with  v  s  w  s  Ef[X{oo)]. 

Expressions  for  the  variance  constants  o;  and  o*  have  been  given  by  Kemeny  and 
Snell  (1960,  1961)  and  Hordijk,  Iglehan  and  Schassberger  (1976);  see  also  G)yim(l984)  and 
Whitt  (1991).  Moreover,  Theorem  4.2  of  Hordijk  et  al.  (1976)  implies  that  ai  <  az 
(Hordijk  et.  al.  ( 1976)  do  not  work  directly  with  the  averages  in  ( 1 )  and  (2)  but  the  subsequences 
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in  which  t  is  rqjlaced  by  the  random  epoch  of  the  n***  visit  to  a  fixed  state,  i.e.,  regeneration 
points.  The  comparison  is  fair  because  both  processes  have  the  same  random  number  of 
transitions.  The  comparison  carries  over  to  (1)  and  (2);  see  Lemma  1  on  p.  210  of  Glynn  and 
Whitt  (1986)  and  §6  of  Glynn  and  Whitt  (1989).  The  comparison  of  (1)  and  (2)  is  fair  because 
both  processes  have  the  same  expected  number  of  transitions.)  Hence,  unlike  Example  1.  here  the 
customer  average  W(t)  is  always  asymptotically  more  efticient.  ■ 

Estimating  the  averages  is  also  important  when  ASTA  is  not  known  to  hold.  We  may  want  to 
use  a  customer  average  as  an  approximation  for  a  time  average,  or  vice  versa.  Then  we  might 
want  to  estimate  the  difference  to  see  if  the  a{^ximation  is  reasonable,  and  possibly  make 
corrections.  For  example,  estimating  time  averages  can  be  computationally  intensive  in  discrete- 
event  simulation.  Time  averages  can  be  estimated  if  we  collect  data  at  every  event  and  use  the 
time  intervals  between  these  events.  However,  tlw  total  number  of  events  may  be  prohibitively 
large.  An  alternative  procedure  is  to  introduce  extra  independent  events  strictly  for  sampling. 
However,  it  is  typically  easier  to  collea  data  only  at  time  points  in  a  subset  of  the  original  events. 
We  might  routinely  use  such  a  potentially  biased  procedure  to  estimate  time  averages,  and  only 
occasionally  do  a  more  complicated  estimation  as  described  here  in  order  to  estimate  the  bias  to 
see  if  it  is  serious  and.  if  so,  make  corrections. 

As  a  theoretical  basis  for  constructing  confidence  intervals  for  estimators  of  the  two  limiting 
averages  and  their  difference,  we  establish  supporting  central  limit  theorems.  In  §1  we  specify 
the  model.  In  $2  we  show  that  the  desired  CLT  follows  from  a  CLT  for  basic  component 
processes.  In  §3  we  introduce  martingale  structure  and  give  a  sufficient  condition  for  the  CLT  in 
§2.  In  §4  we  show  that  the  conditions  in  §3  are  satisfied  in  a  certain  Markov  framework.  For  this 
purpose,  we  apply  the  martingale  CLT  on  p.  339  of  Ethier  and  Kurtz  (1986).  For  additional 
background  on  CLTs  for  general  processes,  see  Chapter  4  of  Billingsley  (1968).  Section  7.7  of 
Durrett  (1991)  and  Section  9.2  of  Lipster  and  Shiryayev  (1989). 
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In  §S  we  show  how  the  martingale  CLT  in  §4  can  be  applied  to  generalized  semi-Markov 
(»t)cesses  (GSMPs).  In  §6  we  consider  the  special  case  of  continuous-time  Markov  chains. 
Finally,  in  §7  we  draw  conclusions. 

This  paper  was  partly  motivated  by  Yao  (1991).  Yao  (1991)  noted  that  the  martingale  CLT 
applies  to  a  martingale  associated  with  PASTA  introduced  by  Wolff  (1982).  However,  this 
martingale  is  not  the  natural  estimator  for  the  difference  of  the  two  limits;  see  Propositions  S  and 
6  below.  The  analyas  here  is  also  similar  in  spirit  to  the  analysis  related  to  indirect  estimation 
via  L  s  A.  W;  see  Law  (197S),  Carson  and  Law  (1980)  and  Glynn  and  Whitt  (1986,  1987,  1988. 
1989).  A  new  paper  which  presents  additional  related  results  concerning  comparisons  of 
estimators  is  Glasserman  (1991). 

1.  The  Model 

Consider  a  continuous-time  stochastic  process  X  m  {X(0  ^  ^  0)  taking  values  in  a  general 
space,  and  a  stochastic  point  process  on  the  interval  [0.  <»).  characterized  by  the  counting  process 
(V  ■  {W(r)  ;  r  ^  0}  or,  equivalently,  the  sequence  of  successive  points  {r„  :  n  ^  1 );  i.e., 
JV(t)  s  sup{n  2  0  :  ^  r},  r  2  0,  where  Tq  =  0  without  there  being  a  O"'  point.  Because  of 

the  many  applications  in  queueing  theory,  we  refer  to  Af  as  the  arrival  process,  even  though  the 
points  need  not  be  interpreted  as  arrivals.  Let  (/(r)  =  f{X{t))  where  /  is  a  real-valued 
measurable  ftinction  on  the  state  space  of  X.  We  assume  that  the  sample  paths  of  (/  are  left 
continuous  with  right  limits,  while  the  sample  paths  of  ^  are  right  continuous  with  left  limits. 

We  are  interested  in  the  time  average  V{t)  in  (1)  and  the  customer  average  W{t)  in  (2),  within 
the  general  framework  here.  Under  appropriate  regularity  conditions. 

V(r)  -♦  V  and  W(r)  w  w.p.  1  as  r  -» oo  (7) 

where  the  limits  v  and  w  satisfy  the  bias  formula 
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w  -  cov[4*(/).  t/‘(/)] 


(8) 


with  U*(t)  =  X*  and  N*  being  the  stationary  and  stationary-increment  versions  of  X 


and  N,  and  m.*  (t)  being  the  conditional  intensity  of  a  point  from  N*  at  t  conditional  on  X*  (r);  see 
Brdmaud  (1989),  Br^maud,  Kannurpatti  aiui  Mazumdar  (1991).  Makowski,  Melamed  and 
Whitt  (1989)  and  Melamed  and  Whitt  (199()a,b).  See  (17)  below  for  an  analog  of  (8)  derived  via 


the  CLTs. 


Moreover,  under  appropriate  regularity  conditions  (see  Sections  2-4),  in  this  general 
framework,  CLTs  (3)  and  (4)  hold  plus  a  joint  CLT  for  ( W(t) ,  V{t)),  i.e., 

f‘'^(V(/)-v  ,  W{t)-  w)  ^  NiO,  I)  as  /  -»  oo  (9) 

where  the  covariance  matrix  £  has  diagonal  elements  (variances)  Oy  and  al,  aiKl  off-diagonal 
elements  covariances  9^.  As  a  consequence, 

r'^(W(r)  -  V(r) -(w- V))  => /V(0,  03)  as  t ->  «  ,  (10) 


where 


=  ai  +  <jI  -  2a»w  .  (11) 

Note  that  the  limiting  averages  v  and  w  appear  in  the  CLTs  in  the  translation  terms. 

In  applications  we  will  thus  want  to  determine  the  variance  constants  Oy ,  o^,  Oy^,  and  al  in 
(3),  (4)  and  (9)-(II)-  In  some  cases  we  will  be  able  to  calculate  the  variance  constants 
analytically,  but  usually  we  must  estimate  them.  For  example,  we  might  estimate  al,  using  batch 
means,  i.e.. 


ai  =  i  (nlW(kt/n)  -  W({k  -  l)i/n)]  -  ■ 


See  Damerdji  (1989a.b).  Glynn  and  Iglehan(1988),  Goldsman  and  Meketon  ( 1986)  and 
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references  cited  there  for  more  information  about  variance  estimators. 

2.  A  General  Central  Limit  Theorem 

We  first  give  general  conditions  for  the  vector  process  [V(/),  W(/)]  and  the  difference 
W{i)  -  ^(/ItoobeyCLTs. 

Proposition  1.  (a)  If  there  exist  constants  v,  w  and  X  such  that  as  t  -*  oo 

^  N  '  U{s)ds  -  vt  ,  U{s)N(ds)  -  Xwt ,  N{t)  -  Xf  »  ^2-  ^3)  •  ^^2) 

then 

—  V  ,  W(t)  —  H’j  ^  (Li,  X~' {L2  -  wLj ))  as  /  »  (13) 

and 

[w(0  -  V(0  -  (w  -  V)  =>  X~‘L2  -  X"*wL3  -  Li  asr-»««.  (14) 

(b)  If,  in  addition,  (L\,  Li,  L^)  is  normally  distributed  with  means  0,  variances  a}  and 
covariances  Oij,  then  the  limits  in  ( 13)  and  ( 14)  are  normally  distributed  with 

Ov  =  o?  ,  ai  =  X~^(oJ  -  2w023  +  w^o^)  and 

02  f  Oi  +  +  X^^W^O?  -  2X"‘Oi2  +  2X"‘h'0,3  -  2X"*>V023  . 

Proof.  Note  that 

WU)  .  i !'  U(sm<is>  -  [^1  ±  j;  uumJs) , 

SO  that 

,''^(*-(0  -  »)  .  U(s)f/m  -  X».]  -  [^]  ^  I' 

^  X~*L2  -  (X"')(L3)(X“‘Xw)  =  X"‘(f-2  “  wi-s)  as  f  -♦  00  , 
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Hiis  convergence  above  jointly  with  the  normalized  version  of  V{t)  follows  from  the 
convergence  in  (12)  Jointly  with 


t  _1_ 

mo  ’  xt 


U{s)mds) 


(A,"*  .  w)  , 


which  holds  by  Theorem  4.4  of  Billingsley  (1968)  and  the  continuous  mapping  theorem. 
Theorem  S.  1  of  Billingsley  ( 1968).  ■ 

Remark  2.1.  The  asymptotic  efficiency  of  the  estimator  [  V(t),  W(0]  in  (13)  can  be  improved 
by  using  a  linear  control  estimator  if  X  is  known;  see  Section  8  of  Glynn  and  Whitt  (1989)  and 
references  cited  there.  Moreover,  if  v  =  w.  then  we  can  further  improve  the  asymjMotic 
efficiency  by  the  same  method. 


3.  Martingale  Structure 


We  now  proceed  to  exploit  martingale  structure;  for  background,  see  Br^maud(1981)  and 
Ethier  and  Kurtz  (1986).  Let  {  :  /  >  0 }  be  the  internal  history  or  filtration  of  (X,  NY,  i.e.,  S’, 

is  the  a-held  generated  by  { [X(s).  fV(s)]  :  0  ^  s  S  r).  Recall  that  we  assumed  that  U{0  is  left 
continuous,  so  that  it  is  predictable.  We  also  assume  that  the  point  process  A^(r)  has  a  stochastic 
intensity  p(r).  (The  CLTs  below  extend  to  point  processes  with  more  general  compensators.) 


We  now  give  a  sufficient  condition  in  this  framework  for  the  condition  of  Proposition  1  to 
hold. 


ProposiUoB  2.  The  CLT  (12)  holds  if  there  exist  constants  v,  w  and  X  such  that 
^  •  JJc/(5)n(5)ds  -  Xwr , 

JJ(/(j)(/V(di)  -  \i(s)ds]  ,  N{t)  -  /Jn(j)<ir]  . 2$)  as  /  -»  <»  . 


(15) 


in  which  case  L\  =  Z\,  Li  =  Z3  +  Z4  and  Z5.  In  addition,  the  asymptotic 
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normality  condition  of  Proposition  lib)  is  satisfied  if{Y\ . YsMs  normally  distributed. 

Proof.  Apply  the  continuous  mapping  theorem  with  addition.  ■ 

Since  U  is  predictable  and  p.  is  the  intensity  of  N,  the  fourth  and  fifth  components  on  the  left 
in  (15)  are  martingales  under  the  usual  moment  conditions.  Hence,  CLTs  can  be  established  for 
these  components  directly,  as  Yao  (1989)  did  for  the  fourth  component  when  p(/)=X.  w.p.l.  (For 
the  CLT  plications,  it  suffices  to  have /bounded  and  E[N(t)^  ]  <  oe  for  all  r.)  The  other  three 
processes  are  centered  cumulative  processes,  which  satisfy  CLTs  under  extra  conditions  such  as 
appropriate  mixing  or  regenerative  structure,  e.g.,  see  Theorem  20  of  Billingsley  ( 1968)  and 
Glynn  and  Whitt  (1987).  We  will  establish  the  joint  convergence  of  all  five  components  in  the 
presence  of  extra  Markov  structure  by  showing  that  the  vector  process  is  asymptotically 
equivalent  to  a  martingale.  This  Markov  structure  plies  to  a  large  class  of  models,  e.g.,  after 
adding  supplementary  variables;  see  Glynn  (1989a). 

Remarks  3.1.  Note  that  nothing  is  gained  by  going  from  Proposition  1  to  Proposition  2  when 
PASTA  holds,  because  then  p(r)  s  \  w.p.l  and  v  s  w;  then  the  second  term  on  the  left  in 
Proposition  2  is  identically  0,  while  the  third  term  is  X  times  the  first. 

3.2.  Note  that  if  the  CLT  (15)  holds,  then  the  limiting  averages  v  and  w  can  be  identified 
ftom  the  asymptotic  behavior  of  [  U(t) .  p(f)),  without  directly  con.idering  N{t).  In  particular,  if 
(15)  holds,  then 

[V(r),  r“'  JJii(s)ds  ,  r'  JJl/(s)p(s)<fal  4  [v.  X.  Xw]  .  (16) 

p 

where  -*  denotes  convergence  in  probability,  so  that 

p 

-»  w  as  /  -»  oo  ,  (17) 

f>(..  ^ 

which  “explains”  (8)  and  complements  Melamed  and  Whitt  (1990b)  and  Brdmaud  (1989). 
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4.  A  Markov  Framework 

In  this  section  we  assume  that  the  basic  processes  X  and  N  are  defined  in  terms  of  an 
underlying  Markov  process  Y,  which  is  assumed  to  have  left-continuous  sample  paths.  As 
regularity  conditions,  let  the  state  space  of  y  be  a  complete  separable  metric  space  and  let  the 
sample  paths  be  continuous  from  the  left  and  have  limits  from  the  rigtu.  Let  the  filtration 
{  :  r  >  0 }  be  the  internal  filtration  of  Y. 

In  particular,  we  assume  that  there  are  measurable  functions  g  and  h  so  that 

XU)  =  g{YU))  .  \iU)  =  HYU))  and  UU)  =  /(g(YU)))  ,  (18) 

where  /is  the  real-valued  function  defined  in  §1  and  is  the  stochastic  intensity  defined  in  §2. 
Cases  in  which  the  intensity  of  Af  is  a  function  of  Y  occur  naturally  when  the  points  of  N  are 
completely  determined  by  jumps  of  Y;  see  Sections  of  Melamed  and  Whitt (1990a),  (3.2)  of 
Melamed  and  Whitt  (1990b)  and  p.  597  of  Serfozo  ( 1989). 

As  a  consequence,  the  three  cumulative  process  in  Proposition  2  can  be  represented  as 

U(s)ds  =  JJil:,(l'(j))ds 

|Jn(5)iir  =  JJk2(y(r))<fa  (19) 

t/(r)tt(s)ds  s  JJ  kj{Y(s))ds 

wiOiki(y)  •  f(g(y)),  k2(y)  =  h(y)mik}(y)  =  /(g(y))fi(y).  We  provide  conditions  under 
which  these  cumulative  processes  are  asymptotically  equivalent  to  martingales.  Let  A  be  the 
infinitesimal  generator  of  Y  and  2i(A)  its  domain.  We  assume  that  there  exist  constants  y,  and 
functions  b,  defined  on  the  state  space  of  Y  and  belonging  to  2)(A)  satisfying  Poisson’s  equation. 
i.e., 

Abi  =  -  kj  +  7,  (20) 

for  i  3  1,2  and  3.  Typically 
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Y,  =  rtjfc,  a  I  ki{x)K{dx)  (21) 

where  ir  is  the  unique  invariant  probability  measure  of  Y.  (Typically  the  Markov  process  will  be 
irreducible  and  positive  recurrent.)  See  Glynn  (1984,  1989b, c)  and  Whitt  (1991)  for  further 
discussion. 

Remark  4.1.  In  this  se.ting,  the  limiting  averages  v,  \  and  Xw  in  (16)  coincide  with  the  steady- 
state  means  Y  i  •  Ya  Yi  obtained  via  ( 19)-(2 1 ).  ■ 

Since  b,  e  3(4), 

M,U)  ■  biino)  -  fj  {Ab,)iY{s))ds  (22) 

is  a  martingale  with  respect  to  the  filtration  of  K;  see  p.  162  of  Ethier  and  Kuitz  ( 1986). 
Propositions.  !f  (20)  holds  and  ^  0  as  t  oo  for  i  =  1,  2,  3,  then 

(ij  ki{Y{s))ds  -  fit  -  Af,(f)]  =**  0  as  r  -»  «  , 

■where  Mi  is  the  martingale  in  (22). 

Proof.  If  (20)  holds,  then 

ki(Y{s))ds  -  Yi^  =  -  {Abi){Y{s))ds 

=  -  biiYiO)  +  biiYit))  -  {Ab,){Y{s))ds  . 

Apply  (22)  and  the  condition.  ■ 

Note  that  if  T  is  assumed  to  be  stationary,  then  one  condition  in  Proposition  3  is  automatically 
satisfied,  i.e., 

d 

r'^b.d'd))  =  r'^b,(P(0))  =*-0  as  /  -»  oo  ,  (23) 

d 

where  =  denotes  equality  in  distribution. 


Hence,  to  establish  the  CLT  (IS),  it  remaias  to  establish  the  condition  of  Proposition  3  and  the 
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conditions  of  a  martingale  CLT.  We  apply  the  martingale  CLT  on  p.  339  of  Ethier  and 
Kurtz  (1986).  For  additional  background  on  CLTs  for  general  processes,  see  Chapter  4  of 
Billingsley  (1968),  Section  7.7  of  Durrett  (1991)  and  Section  9.2  of  Lipster  and  Shiryayev  (1989). 
In  particular,  by  the  convergence-together  theorem  (Theorem  4.1  of  Billingsley  (1968)  or 
Corollary  3.3  on  p.  1 10  of  Ethier  and  Kurtz  (1986)).  it  suffices  to  apply  the  martingale  CLT  to  the 

five-dimensional  process  M(/)  =  . Ms(i)]  having  the  first  three  components  in 

(22),  with  ( 19),  and  the  last  two  components  in  (15),  i.e., 

MiU)  =  f'  U{s)[N(ds)  -  \i{s)ds]  and  MsU)  =  N{i)  -  J'  p(s)ds  .  r  >  0  .  (24) 

0  '0 

Hence,  we  are  in  a  position  to  apply  the  martingale  CLT.  We  state  the  martingale  CLT  for  M  as  a 

condition  in  our  proposition.  Sufficient  conditions  appear  on  p.  340  of  Ethier  and  Kurtz  ( 1986). 

Proposition  4.  If  the  conditions  of  Proposition  3  hold  and  t~  M(/)  =>  Z'  for  M  in  (22)  and 

(24),  then  the  CLT  (IS)  holds  and  Zi  =  Zj,  I  S  i  ^  5. 

Under  the  regularity  conditions  given  on  p.  340  of  Ethier  and  Kurtz  (1986),  /"'^M(f)  is 
asymptotically  normally  distributed  as  /  -»  ee  with  covariance  elements  determined  by  the 
quadratic  variation  processes  associated  with  M.  In  particular,  suppose  that  M  is  square 
integrable,  i.e.,  £(|M(r)|^l  <  «.  Then  we  have 


(A/,.  A/.KO  p  , 

- >  o*  as  ;  oo 

t 


(25) 


{Mi,M,\(t)  p 

-  -»  o,,  as  /  -♦  oo  , 

t  ’ 


(261 


where  [Mi,  Mi\  ■  (A/J  is  the  quadratic  variation  and  [A/,,  Mj]  is  the  cross  variation:  sec 
pp.  67, 79  of  Ethier  and  Kurtz  ( 1986). 

Remarks  42.  Another  convenient  view  of  the  variance  and  covariance  terms  for  1  ^  i  <  i 
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comes  from  considering  a  stationary  process  framework,  which  is  achieved  by  letting  K(0)  have 
the  unique  invariant  probability  measure  n.  Then 

Oy  =  2f“  Cov[ki(Y{0)).kj{Y{i))]dt  - 

e.g.,  see  Theorem  20. 1  of  Billingsley  ( 1968)  and  Section  9.2  of  Lipster  and  Shiryayev  ( 1989). 

43.  The  CLT  on  p.  339  of  Ethier  and  Kurtz  (1986)  is  expressed  in  a  stronger  functional  CLT 
form.  There  a  sequence  of  martingale  processes  {Mn:n>l}  is  defined  with 
Mn(r)  =  M(nr)/V^  r  ^  0.  The  associated  limit  process  is  then  hve-dimensionai  Brownian 
motion.  We  obtain  the  standard  limit  in  referred  to  above  by  applying  the  continuous 
mapping  theorem  with  the  projection  map  at  time  r  =  1.  ■ 

From  (24)-(26),  it  is  easy  to  identify  the  limits  for  4  £  ij  ^  5.  For  this  purpose,  consider  the 
typical  case  in  which 

U{Tt)  W  and  E[U{Ti)‘]  E[W‘]  as  k oo  ,  (27) 


so  that 


Then 


E[U(Tt)]  ->  E[W]  =w  ask-*oo 


(28) 


t 

1  /V(»)  P 

-  t  045  = 

*  *-I 

-  I  ViTt)'  03  =  \E[W'] 

^  km\ 

Moreover,  if/is  an  indicator  function,  then  Udn)'  s  f/( )  and  03  s  O45  «  Xw 


(29) 


Our  next  result  is  a  CLT  for  the  difference  W{t)  -  V(t)  when  PASTA  holds.  Equivalently. 


we  assume  that  the  stochastic  intensity  ji(  i)  is  constant  (so  that  Af  must  be  Poisson). 
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Proposition  5.  In  addition  to  the  conditions  of  Proposition  4,  suppose  that  (7)  holds,  f  is 
bounded,  E[N(t)^  ]  <  <»  and  n(0  =  k  w.p.J  for  all  t.  Then  the  CLT  ( 14)  holds  with  v  =  w  and 

a3  =  A.-‘Var[H']  .  (30) 

Proof.  ASTA  holds  because  (7)  holds. /is  bounded  and  |i(/)  =  X  w.p.l;  e.g.,  see  Theorems  1 
and  2  of  Melamed  and  Whitt  (1990b).  Consequently,  v  =  w,  =  0  and  Zj  =  kZ\\  see 
Remark 3.1.  The  limit  in  (14)  is  X“*L2  -X"'vv'L3  -  L\,  which  coincides  with 
X“'  Z4  -  X“ '  wZs  under  the  extra  assumptions  here.  Thus,  by  (29), 

03  =  X"-05  X-*w2  -  2X-‘h-045  =  X-'  Var(W]  .  ■ 

Example  1  Revisited.  As  a  consequence  of  Proposition  S.  for  the  M/M/1  queue  with  arrival  rate 
equal  to  the  traffic  intensity  p  discussed  at  the  outset, 

al  =  X-‘Var(W]  =  .  OD 

(1  -  pr 

Note  that  in  (31)  is  asymptotically  negligible  compared  to  a;,  and  at  ~  a;  in  (5)  and  (6) 
asp  ^  1.  From (11) and (31), 

>  0  .  .  (32) 

2(1 -P)‘ 

As  a  consequence  of  Proposition  S,  when /is  an  indicator  function,  we  obtain  Yao’s  (1991) 
(XT  for  Af4(t)  plus  a  special  form  for  the  CLT  for  the  difference. 

Propusitloa  €.  In  addition  to  the  conditions  of  Proposition  5.  suppose  that  f  is  an  indicator 
function.  Then 

^  N(0.  Xw)  as  r  00  (33) 

and  (14)  holds  with 


o2  *  X“'w(I  -  w)  . 


(34) 
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Proof.  Apply  (29)  and  (30).  ■ 

As  a  simple  check  of  (33)  and  (34),  let/  be  identically  1.  Then  obviously  W(t)  =  VI/)  =  1 
for  all  t,  so  that  03  =  0. 

Note  that  the  relative  standard  error,  defined  by  is 


which  becomes  large  as  vv  gets  small.  Hence,  the  difference  of  small  probabilities  are  hard  to 
estimate  in  the  sense  of  (3S).  An  example  would  be  small  blocking  probabilities  in  a  loss  model. 

5.  Generalized  Semi-Markov  Processes 

In  applications  to  queues  and  related  stochastic  models,  the  underlying  Markov  process  Y 
often  will  be  associated  with  a  generalized  semi-Markov  process  (GSMP);  see  Glynn  (1989a).  To 
treat  this  case,  suppose  that  the  sample  paths  of  K  are  piecewise  linear  with  jumps  occurring  at  the 
non-differentiable  points  of  the  sample  paths.  Let  {S„  :  n  ^  1 }  be  the  jump  times  of  Y  and  let 
J{t)  count  the  number  of  jumps  of  Y  in  [0. 1].  We  assume  that  all  the  jump  times  of  N  are  also 
jump  times  for  Y,  i.e.,  ( r„  :  n  ^  1 }  is  a  subsequence  of  {5„  ;  n  >  1 }.  Then  the  five  quadratic 
variation  processes  and  the  ten  cross-variation  processes  associated  with  M  are: 

m 

[Mi.MiUt)  =  2  (^(n5t  +  ))  -&.(T(5t)))2  ,  i  =  1.2.  3. 

kmi 

NU) 

=  2  kiiYiTt))-  . 

km\ 

(A/j.A/sKO  =  A^(0  .  (36) 

m 

-  ZibiiYiSfk-))  -  b,{Y{St))]  ■ 

kml 


(1>;(K(S»  +  ))  -  bj{Y{St))\  .  1  ^  i.  ;•  S  3.  i  #  j. 
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S(t) 

[Mi.M4]U)  =  5;  [bimTk  +  ))  -  b,(Y(Tt))]ki{Y(Ti,))  .  1  </ <  3  . 

k^\ 

NU) 

[Mi,Ms]{t)  =  2:  lbi{YiT„  +  ))  ~  b.mTk))]  .  1  <  «■  <  3  . 

k~l 

Nil) 

[M4,Ms](t)  =  2:  *1(^(7*))  . 

What  we  need,  then,  are  weak  laws  of  large  numbers  (WLLN)  as  specified  in  (25)  and  (26)  for 
the  quadratic  and  cross  variation  processes  in  (36).  Of  course,  we  always  have  (29). 

6.  Continuous-Time  Markov  Chains 

As  a  concrete  application  of  §5,  we  now  suppose  that  K  is  an  irreducible  finite-sute 
continuous-time  Markov  chain  (CTMC).  Then  A  is  the  infinitesimal  generator  matrix  and 
Poisson’s  equation  (20)  has  a  solution  for  y,  defined  by  (21)  for  each  1;  see  §4  of  Whitt  (1991). 
Then  all  solutions  of  (20)  are  of  the  form 

bi  =  Zki  +  {nbi)e  ,  (37) 

where  re  is  the  unique  invariant  measure,  e  is  a  column  vector  of  I’s  and  Z  is  the  fundamental 
matrix  with  components 

Zim  ^  jJlPimU) -n„]dt,  (38) 

with  P{t)  being  the  time-dependent  transition  matrix  associated  with  A;  see  (13)  and 
Propositioo  4.2  of  Whitt  (1991).  In  this  setting  the  conditions  for  the  CLT  ate  easy  to  verify. 
Proposition  7.  If  Y  is  an  irreducible  finite-state  CTMC,  then  the  conditions  of  Proposition  4  are 
satisfied  so  that  the  CLT  (15)  holds.  The  limit  has  a  multivariate  normal  distribution  with 
covariance  matrix  elements  given  by  (15),  (26),  (29)  and  (36). 

Given  that  7(5t)  «  /andy(5ff)  a  m. 
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bi(Y(Sk  +  ))  -  biiYiSt))  =  (Zki),  -  {Zk.)„,  =  I(Z,;  -  Z^,)fc,(y)  .  (39) 

J 

From  (36)  and  (39),  we  thus  deduce  familiar  expressions  for  a}  and  o.y  for  1  ^  i,  j  ^  3;  e.g.,  see 
Remark  4.2  above  and  (12)  of  Whitt  (1991).  For  ways  to  calculate  a}  and  a,j,  see 
Proposition  4.2  and  its  Corollary  3  in  Whitt  (1991). 

Proposition  8.  For  the  finite-state  CTMC  and  1  S  i.  ;  <  3. 

of  =  2X  j:ki{l)K,Zunk,(m) 

I  m 

=  21  I,>^ia)i^iZ,„kj(m)  . 

I  m 

Proof.  We  only  derive  the  expression  for  <s};  tl^  derivation  for  Oij  is  similar.  Let  A/  be  the  sum 
of  the  off-diagonal  elements  of  the  /'''  row  of  A.  Since  Ai/Ai  is  the  transition  probability  of  the 
embedded  chain  in  state  /  at  epochs  5^, 

P{Y{St)  =  /)  -»  K/A//J  KjAj  as  *  oo 

J 

and 

-»  2  KjAj  w.p.  1  as  t  oo  . 


Hence,  by  (36)  and  (39), 
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ar  = 


Z  Z  Z 

I  t  nt 


K,A, 

Z  ^jAj 

J 

A, 

[(Zk,),  -  {Zk,),„\‘ 


=  Z  I  i^iA,„[(Zk,),-(Zki)„]" 

I  m 


=  Z  Z  T^iAunKZk,)}  -  2(Zki),(Zki)„  +  (2k, )‘] 

/  m 


=  -2  2  Z  ^/^/<»i(2k,)/(Zk,),„  (explanation  below) 

I  m 

~  Z  Z  Z  Z  ^lAim  Zijki{j)Z^ki{n) 

I  m  J  n 

=  2  Z  Z  ^i(^)^tZijki{j)  (explanation  below)  . 

I  J 

The  fourth  equality  above  holds  because  JtA  =  0  and  Ae  =  0.  The  final  equality  holds 
because 

Z  ~  ~  ^i(0  F 

m  R 

which  holds  by  (20)  and  (37).  ■ 


By  similar  reasoning,  we  can  also  obtain  expressions  for  0,4  and  0,5,  1  S  i  <  3,  if  we  make 
some  addittonal  assumptions.  Suppose  that  the  subsequence  { T* )  corresponds  to  all  transitions 
from  /  to  m  for  a  set  of  pairs  (/,  m).  Then  { K(r*) }  is  a  discrete-time  Markov  chain  (DTMC) 
with  stationary  distribution  say  tc* . 

PropoalUon  9.  //  { T* )  and  it*  are  defined  as  above  and  f“  ’  Af(r)  -*  X.  as  r  -♦  00.  then 
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and 

0.5  =  >.  I  I  I(Zk.)/  -  (».)«.] 

/  m  I 

for  1  ^  S  3. 

7.  Conclusion 

We  have  provided  a  general  framework  for  establishing  central  limit  theorems  for  time  and 
customer  averages.  The  value  of  this  approach  is  demonstrated  by  the  explicit  formula  for  the 
difference  variance  constant  when  PASTA  holds  in  (30).  The  explicit  expressions  for  the 
quadratic  variation  and  cross  variation  processes  in  a  GSMP  framework  in  (36)  should  also  assist 
further  study. 
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In  this  paper  we  establish  a  joint  central  limit  theorem  for  customer  and  time  averages  by 
applying  a  martingale  central  limit  theorem  in  a  Markov  framework.  The  limiting  values  of  the 
two  averages  appear  in  the  translation  terms.  This  central  limit  theorem  helps  to  construa 
confidence  intervals  for  estimators  and  perform  statistical  tests.  It  thus  helps  determine  which 
finite  average  is  a  more  asymptotically  efficient  estimator  of  its  limit  As  a  basis  for  testing  for 
PASTA  (Poisson  arrivals  see  time  averages),  we  determine  the  variance  constant  associated  with 
the  central  limit  theorem  for  the  difference  between  the  two  averages  when  PASTA  holds. 
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